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ABSTRACT 

Roughly 80% of Ultraluminous Infrared Galaxies (ULIRGs) show blue shifted absorption in the res- 
onance lines of neutral sodium, indicating that cool winds arc common in such objects, as shown by 
Rupke et al and by Martin. The neutral sodium (Nal) columns indicated by these absorption lines are 
~ 10^"^ — 3 X lO^^cm"^, while the bolometric luminosity varies by a factor of only four. We show that 
the gas in ULIRG outflows is likely to be in photoionization equilibrium. The very small ULIRG sample 
of Goldader et al. demonstrates that the ratio of ultraviolet flux to far infrared flux varies by a factor 
~ 100 from object to object. While the Goldader sample does not overlap with those of Rupke et al. 
and Martin, we show that such a large variation in ultraviolet flux will produce a similar variation in 
the column of neutral sodium for a fixed mass flux and density. However, if the cold gas is in pressure 
equilibrium with a hot outflow with a mass loss rate similar to the star formation rate, the range of 
ionization state is significantly smaller. Measurements of the UV flux for objects in the Martin and 
Rupke et al. catalogs will definitively determine if photoionization effects are responsible for the wide 
variation seen in the sodium columns. If they are, a determination of the gas density and mass loss 
rate in the cool winds will follow, with attendant improvements in our understanding of wind driving 
mechanisms and of the effects of galaxies on their surroundings. 



Subject headings: galaxies: general - 
galaxies : star burst 

1. INTRODUCTION 
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galaxy, but that this optically thick component is seen to- 
ward ^ 80% of ULIRGS. Given the patchy nature of the 
interstellar medium in most galaxies, the last finding is 
not entirely surprising. 

However, another implication of the observations is very 
surprising. As Figure 1 shows, the mean column averaged 
over both samples (Rupke et al. 2005a; Martin 2006) is 
13.7 in the log, with a standard deviation of 0.66 dex, or 
from lO^'^ cm^^ to 3 x 10^^ cm^^. The actual Na I columns 



Recent observations of the sodium resonance line (NaD) 
doublet in Ultra-luminous Infra-red Galaxies (ULIRGS) 
by Rupke and coworkers Rupke et al. (2002), Rupke 
& Veilleux (2005), Rupke et al. (2005a), Rupke et al. 
(2005b) and by Martin (2005; 2006) show that cool out- 
flows from ULIRGS are common. Martin (2005) finds that 
15 of 18 ULIRGS (83%) possess such flows, while Rupke 
et al. (2005b) find a detection rate of 80 ± 7% for 30 local 
ULIRGs. These results indicate that the outfiow emerges 
in most directions. The absorption troughs, which typi- 
cally extend over ^ 300 km s^ , are not black, which at 
first blush might suggest that the outfiows are optically 
thin in the NaD line. However, the doublet line ratio, 
which is equal to 2:1 in optically thin gas, rarely has that 
ratio in Martin's sample, while Rupke et al. (2005a) find 
optical depths for the weaker of the doublet lines that 
range from 0.06 to 7, with an average of 1.5. This finding 
shows that some sight lines to ULIRG galaxies are opti- 
cally thick in the NaD line, that the optically thick out- 
flow covers only about 25 — 30% of the (optically emitting) 
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vary from 7 x 10 cm~ to 5 x 10 cm" in the sample 
of Rupke et al. (2005a). Since the luminosity, and star 
formation rate of the galaxies in both samples varies by 
only a factor of about 4 (the standard deviation is a fac- 
tor of 2), one might expect the Nh column and hence the 
NnuI column toward the galaxies to vary by a similar fac- 
tor; in simple galactic wind models the mass outflow rate 
increases with increasing star formation rate. 

The Na I column can be used to find the mass loss rate 
from each galaxy, in principle. For example, in a smooth 
fiow, the density is related to the mass loss rate by 



Qr^ li{r)n{r)v{r) , 



(1) 



2 



where /u(r) is the mean molecular weight, n{r) is the num- 
ber density of hydrogen, v{r) is the velocity of the flow at 
a distance r from the center of the galaxy, and Q < Air is 
the global wind cover factor in steradians. For simplicity 
we employ a spherical model, although ULIRGs appear to 
form a substantial fraction of their stars in a kiloparsec (or 
smaller) scale disk Condon et al. (1991). 

The hydrogen column through such a wind is 



Nh = 



Mu 



f 

•I ra 



dr 



(2) 



A very rough estimate is found by taking /i(r) = rUp and 

v(r) = f oe , where the latter is the terminal velocity of the 
outflow, which in the spirit of approximation adopted here 
is taken to be the maximum observed blue shift: 



H 



(3) 



where Tq is the radius at which the absorption line forms. 
Solving for the mass loss rate in terms (as far as possible) 
of observed quantities. 
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Nnu 



■''Na, 



(4) 



where Nnui is the (observed) column of gas phase neutral 
sodium, N^a the total gas phase sodium column, Nh is 
the total hydrogen column, and dj^a is the ratio of total 
to gas phase Na; it is larger than unity since muc;li of the 
Na is locked up in dust grains. Neither of the two terms 
in parentheses nor rfjva are measured, nor is ro well con- 
strained. For solar abundances Nn/N^a = 4.9 x 10^ An- 
ders & Grevesse (1989), while Savage & Sembach (1996) 
find that the gas phase abundance of Na in the cool gas of 
the diffuse cloud toward Q Ophiuci is a factor of djvo = 8.9 
smaller than the solar abundance of Na. There is evidence 
of dust in ULIRG outflows; first, NaD is a resonance line, 
so that it scatters continuum photons rather than destroy- 
ing them. Since the winds are nearly spherical (they are 
seen toward 80% of ULIRGS), any photons removed from 
our line of sight by the NaD transition should appear along 
some other line of sight as redshifted emission; since we see 
no such emission, some other mechanism is removing the 
scattered photons, absorption by dust being the number 
one suspect. Second, there is a correlation seen between 
reddening and NaD absorption toward luminous infrared 
galaxies, again suggesting dust in the outflow Veilleux et 
al. (1995); Heckman et al. (2000). In the rest of this 
paper we will adopt solar metalicities and the Milky Way 
ISM d^a just quoted. Using these values, the typical ob- 
served column Nal ~ 10^^ cm~^ would correspond to a 
hydrogen column Nh ~ 10^"^ cm~^ if somewhat less than 
half the gas-phase sodium were neutral. 

Any local heating, such as might be produced by shocks, 
will lead to lower depletion levels, but by less than the fac- 
tor of ten corresponding to all the grains being destroyed, 
given the evidence for dust cited above. 

The force responsible for expelling the wind is uncertain. 

One candidate is ram pressure from outflowing supernova- 
heated gas, which will entrain cold gas from the interstellar 



medium of the ULIRG. A variant on this is ram pressure 
from hot gas produced by an accreting massive black hole 
(an active galactic nucleus, or AGN). The Compton tem- 
perature of both Seyfert and Quasar nuclei is around lO^'K, 
similar to that from supernova heating. Momentum driv- 
ing from supernovae can be significant, if most supernovae 
explode while surrounded by cold gas. Another candidate 
for momentum driving is radiation pressure on dust grains 
embedded in the cold gas. The radiation could arise either 
from stars or from a central black hole. The radiation pres- 
sure from a black hole can exceed that from a starburst. 
Whether such a centrally driven outflow can produce an 
outflow that will be detected from 80% of the sky is far 
from certain. 

In fact it seems likely that both energy and momentum 
driving operate. ULIRGS are enveloped in x-ray emit- 
ting gas, Arp 220 being an excellent example McDowell 
et al. (2003), Clements et al. (2002); the presence of hot 
gas is a necessary (but not sufficient) condition for en- 
ergy driving of a cool outflow. As we have just noted, 
the winds are most likely dusty. The wind optical depth 
to optical or UV emission is likely substantial, while the 
galaxies themselves are optically thick even to the far in- 
frared. The outflow may well be multiphase; the questions 
to be addressed revolve around which phase carries more 
mass, momentum, energy, and metals into the surrounding 
intergalactic medium. 

To answer such questions, we would like to measure the 
mass loss rate in cold gas. If the Nal columns are taken 
at face value, and the ratio Nn/N^ai is assumed to be 
constant, the mass loss rates in cold gas vary from ULIRG 
to ULIRG by a factor of 50. Is this really true? 

There are niunerous possible reasons for the Nal column 
seen in galactic outflows to vary from object to object. The 
mass loss rate will vary from object to object, so that the 
total column of hydrogen, and hence Na, will vary, even at 
a fixed metallicity. The metallicity of the outfiowing gas 
may vary from object to object. Both these properties are 
likely to depend on the star formation rate (and history); 
the star formation rate, at least, varies over only a small 
range for the samples we are discussing. A more likely 
cause of variation in Nal column is a variation in either 
or both the depletion on to dust grains and the ionization 
state of Na. Neutral sodium is a rather delicate atom, 
easily ionized by ultraviolet radiation. The wind may be 
irradiated either by starlight from the host galaxy, or by 
shocks in the wind. 

In this paper we argue that the outfiow is in photoion- 
ization equilibrium, except possibly for the most UV dim 
systems, and then show that if Martin and Rupke et al.'s 
galaxies have spectral energy distributions like those of 
other nearby ULIRGs that have been observed in the ul- 
traviolet, then it may be that the mass loss rates are more 
nearly equal than the observed NaD columns imply; un- 
less the neutral gas is very dense, n ~ 10^ cm""^, the ion- 
ization fraction, or the ratio of neutral to total sodium, 
(Nal/Na), will vary by a large factor from object to ob- 
ject. This follows from the observed variation of ~ 100 
in the ratio of UV to FIR flux seen in nearby ULIRGs 
(Goldader et al. 2002), and simple photoionization calcu- 
lations, as we show below. 
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2. NEUTRAL SODIUM GAS ABUNDANCES 

The ionization state of the wind depends on the gas den- 
sity, the ionizing flux, and the gas temperature. We start 
by estimating the gas density. We will use canonical val- 
ues for the luminosity LpiR = lO^^L© and terminal wind 
velocity Voo = 300 km s~^. 

To obtain a rough lower limit for the density we will as- 
sume a smooth dust driven wind, which captures the entire 
momentum output of the galaxy along a particular sight 
line. Following Rupkc ct al. (2GG5b) wc let Cq denote the 
fraction of sight lines that pass through a wind and Cf 
denote the local cover factor of Nal gas (along a sight line 
through the wind), so that il/An = CnCf. Rupke et al. 
(2005b) find Cq = 0.8 as noted above, and Cf « 0.4 (their 
table 2) for local ULIRGs, values which we will adopt. In 
that case the mass loss rate is given by 



^FIR 



CnCf. 



(5) 



For our canonical numbers, this yields a mass loss rate 
of Mw ~ 22M0yr~^. The star formation rate needed 

to power the bolomctric luminosity (IQ-'^^Lq) is = 
L/(ec2) K. QOMgyr^i, where e » 8 x 10"'' for a Salpter 
IMF with stellar masses between 0.1 and lOOM©. The 
ratio 1) of mass outflow rate to star formation rate is 

Mu, „ „ c 
-T— = OoG/e 

M. Woo 



0.256 



0.8 



0.4 



8 X 10-4 / V 1 



/ 300 km s" 



6) 



Rupke et al. (2005b), using a sample of 30 low redshift 
ULIRGS, find r? = 0.19l°;^; the error bars do not include 
the large uncertainty in the Nal/Na ratio. 

Estimates of mass loss rates in supernova driven winds 
also tend to yield results comparable to the star formation 
rate, e.g., Martin (1999), although it is unclear if the hot 
gas actually escapes from the galaxy in massive galaxies. 

Combining equations (1) and (5) yields the estimated 
density 

n(r V L ) — ^FiR /^-l 

Using the canonical values quoted above, ^, ■ 
r = 500 pc, we find 



1.38mp, and 



L 



1012J: 



© 



/300km/s 



/lOOkm/s 



500 pc 



We have chosen r = 500 pc since this is (roughly) the size 
of both CO (~ 500 pc, Solomon et al. 1997) and radio 
disks (~ 100 pc, Condon et al 1991) in typical ULIRGs. 
The distance between the disk and the outermost point of 
the outflow will depend on age, and is likely to be much 
larger than 500 pc, but most of the column will be accu- 
mulated at small radii, see equation (3). While we have 
derived this result for a smooth cool outflow, a hot super- 
nova driven outflow with the same mass loss rate will have 
the same density, as can be seen from equation (1). 



If the gas is clumpy, the relevant density will be higher 
than this estimate, for the same mass loss rate. In par- 
ticular, if the cool gas is driven out of the galaxy by ram 
pressure from a hot outflow, the cool gas should be in or 
near pressure balance with the hot wind. In that case, 



iih- 



Th 



(9) 



where the siibscript h denotes supernova heated gas. The 
hot gas has Th ~ lO^K, while the cold gas has T <^ 100, 
as shown by the numerical results described below. If the 
mass loss rate in hot gas is of order the star formation 
rate, the cold gas has n ^ 6 x 10^ cm~^ near the launching 
radius. This is slightly higher than the density of the ISM 
in the inner disks of ULIRGs. Since the ionization state 
of the gas depends on the gas density, it may be possible 
to distinguish between ram pressure driven flows and dust 
driven flows. We will return to this point in §4. 

2.1. Timescales 

With this estimate of the gas density we are in a posi- 
tion to estimate several relevant time scales. The flrst is 
the dynamical time scale, 



Tdyn 



5 X 10 



13 



500 pc 



300 km s" 



s. 



(10) 



In a clumpy flow, the size of the cold clumps should be 
miich less than the size of the galaxy. For example, in 
our own galaxy, observations of cool gas on top of the hot 
outflow GSH 242-034-37 (McClurc-Grigiths et al. 2005) 
indicate that the cool gas has a length scale I of order tens 
of parsecs, while sitting on hot gas 1,500 parsecs above the 
galactic plane. The sound crossing time for such clumps 
is 

ini2/^ ' \/300kms-i\ 

for the surrounding hot gas. These clumps have a col- 
umn density Nh « 3 x 10^^ cm"^, a factor of about ten 
lower than the column through the disk of the Milky Way, 
and a density n = Nh/1 ~ Icm"'^ for the observed size 
I = 10 pc, a density similar to that in the disk as a whole. 

In order to estimate the time scale for the evolution of 
the cold outflowing gas in ULIRGs, we need to estimate 
the length scale of the cold gas. To estimate the size of 
a clump of cold gas driven out of a ULIRG disk by a hot 
outflow, we use the inferred Nh for the Nal carrying out- 
flow, and assume that the initial gas density is equal to 
that in the disk midplane, as is the case in the Milky Way; 
I = Nh/u = 10^1 cm-2/105cm-3 = 10^*^ cm. In Arp 220, 
the star forming disk(s) have r » 100 pc with a disk scale 
height H satisfying H/r ~ 0.1 Sakamoto et al. (1999), so 
l/H P:i 3 X 10~^, a factor of ten or so smaller than ratio of 
0.05 seen in the Milky Way. 

This short length scale has important implications for 
the pressure and hence the density in the Nal bearing gas. 

To work out these implications, wc first show that if the 
hot gas is responsible for pushing the cold gas out of the 
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galaxy, the mass loss rate in the hot gas must exceed that 
of the cold gas by the factor An/V, « 5. We then show 
that the density of the cold gas is n ^ 6 x 10^ cm~^ for 
our cannonical numbers. Gas of this density exposed to 
the radiation from a ULIRG will be neutral, with all the 
gas-phase Na in the form of Nal. 

We assume that the cool gas seen in ULIRG outflows 
comes in roughly spherical clumps of radius /, with col- 
umn Nh ~ 10^^ cm^^. These clumps are pushed out of 
the galaxy by ram pressure from the hot outflow, which 
exerts a force Murray et al. (2005) 



(12) 



where Vh is the velocity of the hot outflow. This force must 
exceed the force of gravity. 



r 



(13) 



where Vc is the circular velocity of the galactic disk; we 
note that Vc <, Vh- Introducing the cover factor flh of 
the hot gas (which should be close to An since the sound 
speed of the hot gas is comparable to the bulk velocity of 
the gas), we find that the mass loss rate of the hot gas 
must exceed 



Vh 



(14) 



if the hot gas is to accelerate the cold clouds out of the 
galactic potential. 

The mass loss rate in hot gas must exeed that in cold 

gas by the ratio of the cover factors; we expect that this 
is roughly the inverse of C/Co, or about a factor of five. 
For ram pressure driving to be important, the density of 
the hot gas should be that given by eqn. (8). 

The hot gas will travel accross any embedded cold clump 
in a time 



Tft « 3 X 10** 



1016 



cm 



300 km s"^ 



Ch 



(15) 



much shorter than the dynamical time. 

In a supernova driven fiow shocks may be driven into 
the cold clouds. Let us suppose that the mass loss rate 
in hot gas is equal to the star formation rate (and thus a 
factor ~ 5 times the observed mass loss rate of the cold 
gas in ULIRGs). The pressure in the cold clouds is 

= (to^) (tWk) (16) 

while the pressure in the hot gas is 

^''-«><10"'(6^)(iIk)^^-^-"^' (1^) 

This over-pressured (from the point of view of the cold 
gas) hot gas moves at its sound speed, which is a factor of 
100 or more higher than the thermal sound speed in the 
cold gas. A cold gas cloud that is exposed to the hot gas 



will experience a six-fold increase in external pressure in 
a time of order one percent of the cloud dynamical time, 
i.e., it sees what is effectively an instantaneous jump in its 
external pressure. The result will be a shock of Mach num- 
ber ~ \/6 propagating into the cloud. If the initial density 
of the cold gas is smaller, the shock will be stronger. 



These shocks will cross the cloud in a time 



5 X 10^° 



I 



1016 cmy/ V 2km s 



(18) 



where Vg — y^Ph/rripn is the shock velocity. Since this is 
much less than the outflow time given by equation (10), 
the cold clumps should be in pressure equilibrium with 
the hot gas, if the hot gas is dynamically important; this 
statement is true even for / a factor of 10 or 100 larger 
than our flducial value of lO^^ cm. Indeed, in our galaxy, 
McClure-Griffiths et al. (2005) conclude that the clumps 
they see are in pressure equilibrium with the surrounding 
hot gas. 

We have assumed that the shock is isothermal. In fact, 
the cooling time for the postshock gas is: 



Tcool ~ 1-5 X 10^ 



7 / T, \ /3 X 10-26 erg-icm^s-i 



103 Ky V ^ 

3 X lO^cm-^ 



s, (19) 



where we have scaled to a post-shock temperature Tg ap- 
propriate for a 2 km s"^ shock. The numerical value for 
A is from Dalgarno & McCray (1972). The cooling length 
TcooiVs ~ 4 X IQi^ cm, much shorter than any other length 
in the problem, so the shock is effectively isothermal. The 
cooling column, TcooiVsn « lO^* cm-^, consistent with the 
general result given by McKee & Hollenbach (1980). In 
fact this result is more general. Any shock propagating 
through dense gas will cool on a time scale short com- 
pared to the outflow time at r = 500 pc. The maximum 
cooling column is McKee & Hollenbach (1980) 



N, 



cool 



10i9-l0^"cm 



20 . 



(20) 



for Vs < 20 km s ^ , and less for higher shock velocities. 

The minimum value for A « 10^^^ erg--"- cm^ s"! Dal- 
garno & McCray (1972) for T ^ 30K, so as long as 
n > lOcm-3 the gas will cool to 30K in less than the 
outflow time. 

This result shows that, in the event that the initially 
cool, dense component of the ISM is shocked and hence 
heated in the wind, it will quickly cool again. Since the 
ratio of temperatures of the hot and cold gas exceeds 
10^, while, as we have just shown, the two are in pres- 
sure equilibrium, the density of the cold gas must exceed 
n — 6 X 10'^ cm^'^, as we stated above. 

If the mass-loss rate in hot gas is less than the star 
formation rate, so that some mechanism other than ram 
pressure is responsible for expelling the cold gas from the 
galaxy, Ph will be similar to or smaller than the pressure 
in the cold clouds; the only shocks present will be those 
driven into the clouds by the ram pressure of the hot gas 



5 



as it encounters cold gas. The cold gas will not be pressure 
confined by the hot gas, so the cold gas will expand at its 
sound speed as it flows out. If the initial length scale is 
10^^ cm, and the mass loss rate is that given in equation 
(6), the cold gas will exand to a cover factor C/ of order 
unity by the time it is 500 pc from its launching radius. 

The sodium ions in the post-shock gas will also quickly 
recombinc once the gas cools, as we show next. The re- 
combination time for sodium in cool T ~ lOOK — lOOOK 
gas is 



OiNane 



0.06 cm-3 



Tie. 



s. 



(21) 



In this expression, aNa is the total recombination coeffi- 
cient of Nail to Nal. The radiative recombination coeffi- 
cient aNa{T = lOOK) w 10""cm-''s"i Verner & Ferland 
(1996), which we have assumed dominates the recombi- 
nation. We have scaled to the minimum electron density, 
for gas of any density, found in our numerical photoioniza- 
tion calculations using CLOUDY, described below. This 
minimum electron density occurs in gas with low density 
(n = 6 cm-3 and low UV/FIR flux ratios, ~ 10"^; the 
electron density for SEDs with higher UV/FIR flux den- 
sity ratios is much larger, leading to shorter recombination 
times. 

For warm ionized gas, T k, lO^K, rig w n, while 
a{T — 10'') « 2 X 10^^'^ cm"^ s~^. Thus, for any tempera- 
ture between lOOK and lO^'K the recombination time is far 
shorter than the dynamical time. We conclude that shocks 
are not an efficient mechanism for maintaining sodium 
atoms in a singly ionized state 

As the gas flows away from the disk plane the density 
and pressure in the hot component will decrease as 1/r^, 
so the density in pressure confined cold clumps will also de- 
crease, lengthening Tcool- However the cooling time will re- 
main short compared to the outflow time, which increases 
as r, out to megaparsec distances. 

Even for a smooth flow not in pressure equilibrium with 
a massive hot outflow, the density n « 6, so the cool- 
ing time is TcooI ~ 8 x lO"'^ s, still short compared to the 
outflow time. However, if shocks occur at large enough 
distances, the cooling time may exceed the outflow time. 
Nominally this would occur at distances of order 5 kpc, 
but this neglects the compression of the gas in the shock. 

The photoionization time is 



A'Kr'^{hv) 

dNalLuV 

r 



10° 



500 pc 



-'FIR 



-FIR 



s, (22) 



where a^ai ~ 10~ cm Verner ct al. (1996) is the pho- 
toionization cross section and Lf/y ~ i^L]^^ is the luminos- 
ity at the sodium edge. The photoionization time ranges 
between « 10^ and 3 x 10" s for observed UV/FIR flux 
density ratios (see Figure 2 for the flux density ratios). 

We have seen that both the recombination time and the 

photoionization time arc shorter than the dynamical time, 
at least for temperatures below lO^K, while the post-shock 



cooling time is shorter still. We conclude that the Nal 
laden gas is in photoionization equilibrium. 

2.2. Photoionization equilibrium 

Consider a dilute gas of number density n irradiated 
by the light from a ULIRG. We saw above that both the 
photoionization and recombination time scales arc shorter 
than the dynamical time scale, so the ionization state of 
sodium is controlled by the flux near 5.139 eV, the ioniza- 
tion potential of Sodium I, and by the flux near the Lyman 
edge at 13.6 eV. 

The Lyman edge flux controls whether the gas as a whole 

is neutral. In highly ionized gas, the neutral fraction of hy- 
drogen is given by Ostcrbrock (1989) 



n ghc 



OtH 



(23) 



where ~ 10~ cm^ s~ is the recombination coefficient 

for hydrogen, an ~ 6 x 10~^* cm^ is the photoionization 
cross section, c is the speed of light, and 



(24) 



is known as the ionization parameter. The rate of emission 
of ionizing photons, Lo/{hu) is defined by 



f 



dv 



hu ' 



(25) 



where is the frequency at the Lyman edge, and L^, is 
the flux density (luminosity per Hertz, erg s~^ Hz~^), at 
frequency v . 

The Lyman edge luminosity Lq ~ uoL,^^ is not known 
for ULRIGs, but it is likely to be much less than the flux at 
longer wavelengths, for two reasons. First, most ULIRGs 
are powered primarily by stars, which have large intrin- 
sic Lyman edges. ULIRGs powered by AGN may have 
larger Lyman luminosities; if so they should show less ab- 
sorption in Nal. Second, the dusty gas both in the galaxy 
proper and in the outflows will have a substantial (perhaps 
dominant) fraction of neutral hydrogen, which will absorb 
radiation beyond the Lyman edge. With our canonical 
values for mass loss rate (and hence density), and for ex- 
pected Lyman edge fluxes {Lq <; 10'*^ erg s^'^, correspond- 
ing to smooth extrapolations from UV flux density ratios 
L°/L^^^ w 10~^, see Figure 2), the neutral fraction of 
hydrogen is at least 10""*; the column through the ULIRG 
wind is of order lO^^cm"^, so the wind has an optical 
depth of order unity at the edge. 

Galaxies with no intrinsic Lyman edge and large UV to 
FIR flux density ratios 



-UV , 



/Lr« > 10- 



(26) 



will then have ionized winds, independent of any possible 
heating due to shocks or supernova driving. Such winds 
will not contain an observable column of neutral sodium. 

In contrast, galaxies that have no intrinsic Lyman edge 
but which have flux density ratios less than 10~^ will have 
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winds that are optically thick at the edge, as will almost 
any galaxy with an intrinsic Lyman edge. In such galaxies 
the wind will consist primarily of (nearly) neutral hydro- 
gen, assuming that it is not strongly heated by supernovae 
or an AGN. Such galaxies may or may not have substantial 
columns of neutral sodium, as we now show. 

In gas in which hydrogen is nearly neutral, the neutral 
fraction of sodium is given roughly by 



nNall CLNalC U 



U 



where 



UNal = 



(27) 
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is the ratio of the number density of Nal ionizing photons 
to the number density of hydrogen. The subscript ed refers 
to the sodium photoionization edge, i.e.. X^d ~ 0.2/im, and 
the quantity in angular brackets is the average energy of 
ionizing photons. The numerator in this definition is, by 
analogy with equation (25) 
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Combining these expressions, we estimate the neutral 
fraction of Na as 
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The factor e^'^ accounts for the possible presence of dust. 
We remarked above that the lack of NaD emission lines 
together with the correlation of NaD absorption with red- 
dening suggests that the outflows are dusty, but, as also 
noted above, the outflows do not cover the galaxy. The 
light that we see is thus a mixture of starlight directly from 
the galaxy combined with light filtered through the dusty 
wind. The sodium ions, on the other hand, are embed- 
ded in the wind, so they see the starlight after it has been 
filtered through the dust in the wind interior to their posi- 
tion. To account for this, the flux should be reduced 
from that observed by a factor ~ e"'^'', where is some 
appropriate fraction of the dust optical depth through the 
wind. The hydrogen column through the wind is uncertain 
(since the conversion from Na to H columns is uncertain), 
but if the mass loss rate is of order the star formation rate 
the hydrogen column is of order 5 x 10^^ cm~^ (see equa- 
tion 3). Assuming a galactic gas to dust ratio, the dust 
optical depth would be ^ 2 in the UV. Numerically, 
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where we have taken = 1. In terms of UnuI, 
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This expression is valid for [/jva/ ^ 10"^ . For values of 
the ionization parameter less than this, the neutral sodium 
fraction will saturate (at a value of order 0.5, see the nu- 
merical results below). 



We will see that if the cool outflow is driven by a hot 
wind with a mass loss rate comparable to the star for- 
mation rate, and the cool component has a similar mass 
loss rate, UngI ranges from 10"^ to 10~^; for a dust 
driven outflow with the same mass loss rate the range is 

10-3 ^ UNal < 10. 

The object to object variation in the UV flux is most 
naturally explained by variations in the dust optical depth 
along our line of sight; thus it may be appropriate to set 
Td to some smaller value in objects with (relatively) large 
UV fluxes. On the other hand it is possible that even in 
objects with large UV fluxes the hydrogen column through 
the wind is still of order 10^^ cm~^, but that some of the 
star formation is unobscured; even in the object with the 
highest UV flux in the Goldader et al. (2002) sample, the 
UV flux is less than a tenth of the bolometric flux. The 
mass loss rate may or may not be affected by this bolo- 
metrically small UV flux. If it is not, then w 1 would be 
appropriate even for those objects with large UV fluxes. 

We were not able to find photometry in the UV bands 
for the IRAS sources in Martin (2005) or in Rupke et al. 
(2005a). As a proxy, we used the SEDs for for the luminous 
infrared galaxies observed by Goldader et al. (2002). The 
data, obtained from NED, are shown in Figure 2. The flux 
density near the Na I edge (about 0.2 microns) is smaller, 
relative to the bolometric luminosity, in Arp 220 than in 
IRAS 15250-1-3609. The difference is roughly a factor of 90; 
The relative Na I edge flux density in VV 114 is another 
factor of 4 larger than that in IRAS 15250+3609. The 
other objects in Goldader et al. (2002), IRAS 08572-^3915, 
Mrk 273, IC 883, and IRAS 19245-7245, have flux ratios 
intermediate between IRAS 15250-^3609 and Arp 220. We 
note that Mrk 273 is an AGN, with a luminosity between 
a Seyfert and a quasar. Despite this fact, it has a modest 
UV flux. 

Since the Na edge flux densities of our proxy galaxies dif- 
fer by a factor of ^ 90 (for the same FIR flux) we expect a 
similar difference in the fraction of neutral Na in their out- 
flows. For VV 114, with L^^ /L^^" » lO"'', equation (30) 
predicts a neutral fraction of ~ 2 x 10~^(n/l cvD~^)e'^'^; for 
Arp 220, with L^^ /L^^^ w 2 x 10"'^, the predicted neu- 
tral fraction is 0.01(n/l cm~^)e'^'*. 



3. NUMERICAL CALCULATIONS 

We have carried out numerical calculations of the pho- 
toionization state of gas irradiated by sample ULIRG 
spectral energy distributions (SED). These calculations 
were performed with version 06.02.09a of Cloudy, last de- 
scribed by Ferland et al. (1998). For the sake of simplic- 
ity we modeled the wind using a constant density slab, 
with six different densities (6, 60, 600, 6000, 6 x 10** and 
6 X lO^cm"^), with the lowest and highest densities cor- 
responding to a smooth cold outflow and to a clumpy cool 
medium in pressure equilibrium in a hot outflow with mass 
loss rate equal to the star formation rate. The intermedi- 
ate densities model clumpy cold outflows not driven by hot 
gas. Dust grains were included in the model; we used the 
ISM option of cloudy, which depletes the gas abundances 
so that the total metallicity, gas plus dust, is solar. The 
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dust and gas are mixed in the slab, which hes between the 
continuum source and the observer. 

For each density the cold gas is illuminated by various 
SEDs, all having a total luminosity of 10^"^ Lq, appropri- 
ate for ULIRGs. Three of the SEDs employed were those 
of VV 114, IRAS 15250+3609, and Arp 220. In addi- 
tion we calculated a series of models based on the SED 
of Arp 220. In these models we adjusted the flux density 
at 5.1cV (w 2400 A), varying the FIR to UV ratio from 
L]^^ /L^^^ = 10~^ to 10~^. In doing so we fit a power 
law from the last IR data point to 2400 A. From 2400 to 
1450 A we assumed that i^Li, was flat, i.e., Li, ~ i/ ^ as 
in the 2kpc data of Goldader et al. (2002) for Arp 220. 
From 1450 A to the Lyman edge we simply continued the 
spectrum with the same slope, ^ u^^. 

We note that the UV slope for many of the galaxies 
in Goldader et al. (2002) varies with the aperture size, 
and that the total flux scales roughly with the square of 
the aperture size, the latter indicating that the UV flux 
emerges from a much larger region than the radio or CO 
emission. If the same holds true for the ULIRGs in general, 
the UV flux will not decrease as rapidly with distance from 
the galactic center as will the density in a smooth outflow. 
The result will be an increasing ionization parameter with 
increasing radius; if the cool gas is clumpy, and not pres- 
sure confined, then the density will not decrease as 
and the ionization parameter will not increase as rapidly 
as in a smooth wind. We have not attempted to model 
such effects. 

For models with no intrinsic Lyman edge, the spectrum 
from 13.6 eV to the lone x-ray point at 1 KeV was a simple 
power law. In models with an intrinsic Lyman edge, the 
spectrum from the edge to the x-ray point was also a power 
law, but with a discontinuity at the energy corresponding 
to the Lyman edge. 

We explored other SEDs. We found that the details of 

the continuum shortward of 2400 A do not matter as long 
as the Lyman edge flux is small enough; if the fraction of 
neutral hydrogen is larger than about 10"** {U < 10~^) 
at the face of the slab, the wind is optically thick at the 
Lyman edge and the Nal fraction is controlled by the flux 
at the Nal edge. If this condition on the neutral hydrogen 
fraction is not met, there is essentially no neutral sodium 
unless an unphysically large column {Nh > lO^^cm"^), 
sufficient to produce a Lyman edge, is used. 

Calculations with CLOUDY verify the simple argu- 
ments made in the previous section; the neutral fraction 
of Na can vary by factors of 100 or more due to variations 
in the SED of the host galaxy. The neutral fraction varies 
from a maximum near 1 for an SED with a UV to FIR flux 
density ratio of 10~^ to a low of 10""^ for the spectrum of 
VV 114 (assuming a substantial Lyman break; if there is 
no break, the neutral fraction of Na can be much smaller). 

Figure 3 shows the input spectrum of a CLOUDY run, 
the spectrum transmitted through the outflow, and the 
sum of the transmitted and emitted spectrum. The UV to 
FIR flux density ratio is 10"'' at the face of the slab, and 
there is no Lyman edge. After passing through the slab, 
the spectrum exhibits a strong Lyman edge, and strong 
resonance lines, including Lyman alpha. The effect of dust 



attenuation can be also seen clearly; this attenuation en- 
hances the fraction of Na I over what is seen in calculations 
with no dust. The dashed lines show a similar calculation 
in which an intrinsic Lyman edge appears in the input 
spectrum (which has been shifted downward to distinguish 
the two cases). Note the lack of prominent emission lines, 
since there are few ionizing photons to power such emis- 
sion. 

Figure 4 shows the ratio NtMai/NMa as a function of 
the UV (evaluated at 2000 A) to FIR (100 micron) flux 
density ratio for a range of densities. Wc have assumed 
that the flux shortward of the Lyman edge is a factor of 
10 lower than on the long wavelength side. Prom the last 
UV data point to the edge we assumed a power law of the 
form Li, ~ from just shortward of the edge to the 

first x-ray data point we assumed a simple power law. 

The open points correspond to the artificial SEDs based 
on the SED of Arp 220, assuming a Lyman edge. The 
density is held constant at the value n = 6 cm^'^ (open 
triangles), n = 60 cm^'' (squares), 600 cm^'' (pentagons), 
6000 cm^"^ (hexagons), 6 x 10'' cm^"* (circles), and 6 x 
10^ cm~^ (filled squares); the open triangles corresponds 
to a completely smooth wind with = 22M0/yr, the 
filled squares to pressure confined cold clouds in hot gas 
with a mass loss rate equal to the star formation rate. 

The crosses correspond to the artificial Arp 220 SEDs, 
assuming no Lyman edge, and a smooth wind. The lower 
neutral sodium column in the no-edge case compared to 
the intrinsic-edge case arises as follows: Lyman continuum 
photons penetrating the inner face of the slab keep hy- 
drogen ionized until the optical depth at the Lyman edge 
reaches unity; beyond that depth hydrogen recombines, 
and the neutral fraction of sodium is essentially the same 
as in the case of an SED with an intrinsic Lyman edge. 
When the "Stromgren depth" is larger than the column in 
the wind, as when the UV/FIR ratio is 10~^, the neutral 
fraction of sodium is tiny (~ 10~^), so the corresponding 
cross in Figure 4 is off the bottom of the plot. 

The dashed line in Figure 4 gives the prediction of equa- 
tion (30), which agrees well with the more complete calcu- 
lation for UV/FIR ratios larger than 10"'', in the case 
of an intrinsic Lyman edge. 

The three filled triangles correspond to the measured 
values of Nal for Arp 220, Mrk 273 (data for both taken 
from Heckman et al. 2000), and IRAS 23365-1-3604 (Nal 
column from Martin 2005), as labeled in the figure. In 
each case we assumed that the mass loss rate from the 
host galaxy was given by equation (6), then use the ob- 
served Nal column to find the ionization fraction. The 
error bars reflect only the variation in 77 found by Rupke 
et al (2005b), i.e., only the apparent variation in the mass 
loss rate for objects of the same luminosity. To estimate 
the (miknown) flux at the Nal edge for IRAS 23365+3604 
we extrapolated from the nuclear (2.5 kpc) U-band mea- 
surement of Surace & Sanders (2002), assuming a power 
law extending from the last two optical data points. This 
extrapolation is highly uncertain, as a glance at Figure 2 
shows. 

With only three data points, those corresponding to Arp 
220, Mrk 273, and IRAS 23365, the interpretation of Fig. 
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4 is uncertain. In particular, there is no clear trend of Na 
column with UV flux. The three data points correspond 
to densities in the range 60 — 6000 cm~^, suggesting that 
the flow is slightly clumpy, but not in pressure equilibrium 
with a hot flow of the same total mass loss rate. Clearly 
more measurements of near UV fluxes for galaxies in either 
or both the Rupke et al. and Martin samples are called 
for. 

Figure 5 is similar to Figure 4, except that the ion frac- 
tion plotted is that for Magnesium II over the total Mag- 
nesium abundance. The bulk of the magnesium is singly 
ionized for n < 6 x 10^ cm~^ ; even for higher densities most 
of the atoms are singly ionized for UV to FIR flux density 
ratios larger than 10~^. The near-UV resonance transition 
of Mgll should be a good tracer of neutral outflows; if the 
Mgll column is seen to anti-correlate with UV flux, while 
the Nal column does not vary with UV, it would suggest 
that the cool gas has a high density and is therefore prob- 
ably driven out by hot gas. If the Mgll column does not 
vary with UV flux, while the Nal column does, then the 
gas is most likely low density, and not driven out of the 
galaxy by hot gas. 

4. DISCUSSION 

We have identified two plausible sources of variation in 

the ionization fraction of Nal toward ULIRGs, variations 
in the (Nal) ionizing radiation, and global object-to-object 
variations in gas density. The limited observations of UV 
fluxes from ULIRGs (essentially the objects in Goldader 
et al. 2002) suggest that the former source of variation 
may well be happening. 

Variations in the globally averaged density of outflowing 
gas translate into variations in filling factor of dense gas 
if we assume that does not vary much from object to 
object. The filling factor might vary, if, for example, the 
cold gas were in pressure equilibrium with hot gas, but the 
mass loss rate of hot gas varied from ULIRG to ULIRG. 
Another possible source of variation in filling factor might 
be global differences in the outflow, leading to stronger or 
weaker shocks and hence higher or lower post-shock den- 
sities. However, because the sound crossing times of the 
clumps arc smaller than the outflow time, this source of 
variation does not appear to us to be very promising — the 
clumps would simply re-expand until they reach pressure 
equilibrium. 

Another possible source of variation is variable deple- 
tion onto dust grains. If high velocity shocks pass through 

the gas they will destroy the dust grains, so that no Na is 
depleted onto dust. As we have indicated above, there is 
observational evidence that dust survives in the outflows, 
so perhaps any shocks that have passed through the gas 
are not strong. On the other hand, if the hot gas has a 
high pressure, the density of the cold gas may well be high 
enough to regrow grains, or to cause further depletion of 
Na onto any grains that survive. We have not attempted 
to work out the physics of these processes. 

The referee kindly pointed out that observations of NGC 
6240 Veilleux et al. (2003) show emission lines consistent 
with shock ionization. As we noted above, the column of 
the shock-illuminated and -ionized gas is small compared 



to the total column in the outflow. The emission lines may 
be dominated by shock ionization, but it does not follow 
that the bulk of the gas in the wind is. A clear way to 
tell the difference between photoionized and shock-ionized 
flows is to look for a systematic variation of Nal columns 
with UV flux. If such a correlation is found, the gas is pho- 
toionized. If there is no such correlation, then it is shock 
ionized. (The scaling of UV flux with aperture size found 
by Goldader et al. 2002 shows that in their galaxies at 
least, the UV emission is not dominated by that produced 
in shocks in a wind). 

5. CONCLUSIONS 

The mass loss rates of ultraluminous infrared galaxies, 
as determined by observations of the NaD lines, are be- 
lieved to be substantial — of order one fifth of the star for- 
mation rate. This estimate assumes that the ionization 
fraction of sodium is ~ 0.1. The gross physical proper- 
ties of ULIRGs, including their stellar masses, gas masses, 
sizes, and star formation rates vary by factors of only a 
few, so at first glance the assumption of fixed Nal ioniza- 
tion fraction appears reasonable. However, one property 
of ULIRGS that appears, based only on a few galaxies, to 
vary by a large amount is the ultraviolet flux. We have 
shown that cool ULIRG winds are likely to be in pho- 
toionization equilibrium. We used a simple argument to 
show that the large objcct-to-objcct UV flux variations 
seen in the small number of galaxies that have been ob- 
served in the ultraviolet would be enough to produce sim- 
ilarly large variations in Nal column density, if the density 
of the cool gas is n < 10'' (500 pc/r)^ cm~^. More complete 
photoionization calculations using CLOUDY produce the 
same magnitude of variation in Nal/Na with variation in 
UV flux. They also show, for these moderate densities, 
that the Mgll column should not vary appreciably with 
changes in UV flux. 

If, however, the cold gas is in pressure equilibrium with a 

hot wind having a mass loss rate similar to the star forma- 
tion rate of a typical ULIRG, the cool gas will have a den- 
sity n ^ 10^(500 pc/r)^) cm~^. In that case the ionization 

parameter is UnuI < 10~^, and there will be little varia- 
tion in Nal column from object to object. As a check, the 
(so far unobserved) Mgll columns arc predicted to show 
variations that are anti-correlated with the UV/FIR ratio. 

We conclude that the observed factor of ^ 100 variation 
in the UV/FIR flux density ratios in ULIRGs with simi- 
lar bolometric luminosities is consistent with the factor of 
~ 40 variations in the column density of Na I outflows seen 
by Martin (2006) and Rupke et al. (2005a); the different 
UV fluxes may well be the caiisc of the variation seen in the 
Nal columns, if the cool gas has a density n <J 10"'' cm^^. A 
tight correlation between Na I column and UV flux would 
indicate that the cool gas is not accelerated by a hot wind. 

The origin of the object-to-object variation in the 

UV/FIR ratio is not clear. The known correlation between 
reddening and Nal column in lower luminosity galaxies 
suggests that dust in the outflow may play a role. However, 
the Nal optical depth toward luminous infrared galaxies in 
the sample of Rupke et al. (2005b) actually exceeds that 
toward their (more luminous) ULIRGs. It is unlikely that 
the UV/FIR ratio is smaller in the less luminous objects. 
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Thus dust absorption in the outflow is probably not the 
sole reason for the low UV fluxes in ULIRGS. 

Whether the galaxies studied by Martin (2005) and by 
Rupke et al. (2005a, 2005b) actually have large object- 
to-object variations in their near-UV fluxes is not known. 
We note that the results of Goldadcr ct al. (2002) show 
that near IR or even optical flux measurements are not a 
good proxy for the UV flux; see the data for W114, IRAS 
15250+3609, and Mrk 273 in Figure 2. This demonstrates 
that optical observations will not be sufficient to answer 
this question. Observations by GALEX or by the ACS 
on the Hubble Space Telescope could answer this question 



definitively. 
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Fig. 1. — Neutral sodium (Nal) column as a function of galactic infrared luminosity. Data taken from Rupke et al. (2005a) (triangles) and 
Martin (2006) (squares). The Martin results are for sightlines toward the respective galactic centers. Note the large spread in the column for 
galaxies of approximately the same infrared luminosity. 
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Fig. 2. — Spectral energy distributions of various ULIRGs, normized to unity at lOO/^im; data from the NASA/IPAC Extragalactic Database 
(NED), supplimented by Goldader et al. (2002) (we plot the 2kpc aperture data) and Martin (2005). The solid symbols represent objects 
from Martin (2005); triangles, IRAS 10565+2448; squares, IRAS 17208-0014; pentagons, IRAS 19297-0406; hexagons, IRAS 23365+3604; 
septagons, IRAS 00153+5454. The open symbols are from objects with HST UV observations by Goldader et al. (2002); triangles, VV 114; 
squares, IRAS 08572+3915; pentagons, IC 883; hexagons, Mrk 273; septagons, IRAS 15250+3609; octagons, Arp 220; and circles, IRAS 
19254-7245. The UV/FIR flux density ratio varies by a factor of ~ 100 between Arp 220 and IRAS 15250+3609, and by ~ 300 from Arp 220 
to VV 114. The vertical solid line is at the wavelength of the Nal ionization edge. 
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Fig. 3. — Input spectra of a CLOUDY model with a UV to FIR flux density ratio of 10^'' and no Lyman break (top solid lines) and with 
a Lyman break (top dashed curve); the latter has been displaced downward by a factor of one hundred to clarify the presentation. In both 
cases the middle curve shows the transmitted plus emitted spectrum, while the bottom curve shows the transmitted spectrum alone. 
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Fig. 4. — Neutral sodium fraction, Npfal /^Na, as a function of L^^/L^^^. The open polygons arc calculated using the SED of Arp 
220, but with the near UV data points altered to give various ratios of Luy / Lpijj. The (assumed constant) density of the absorbing gas is 
n = 6 cm~^ (open triangles), corresponding to a smooth wind with = 22Mq/ yr at a distance ro = 500 jjc from the source, n = 60 cm~^ 



(circles). The filled squares give the ionization 



(open squares), n = 600 cm~^ (pentagons), n = 6000 cm^"^ (hexagons), and n = 6 X 10* 

fraction for n = 6 X 10® cm~^, corresponding to cold gas in pressure equilibrium with a hot outflow with a mass loss rate equal to the star 
formation rate. The crosses correspond to n = 6 cm"-^ and various flux density ratios, using SEDs with no Lyman edge. The dashed line is 
the prediction of equations (3) and (31) for n = 6cm~^ for an SED with a Lyman edge. The solid points correspond to the measured Nal 
column for Arp 220, Mrk 274, and IRAS 23365+3604, under the assumption that Mw = JjM©, with = 0.25, and ro = 500 pc. 
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Fig. 5. — The ratio of Magnesium II to total Magnesium, NmqII /J^Mg, as a function of / L^^^. The SED's and symbols are as in 
Fig. 4. The bulk of the magnesium is singly ionized for all but the highest densities plotted. The column of Mgll toward a ULIRG would be 
a decreasing function of Luy /Lpiji, particularly for pressure confined high density clouds. This behavior is the opposite of that predicted 
for the Nal/Na ratio for similarly high density gas. 



